rm(list=ls())

# load necessary libraries
library(plyr)
library(boot)
library(parallel)
library(doParallel)
registerDoParallel(cores=4)
library(ggplot2)

raw.data <- read.csv(file="woc2-analysis-final.csv", stringsAsFactors=F)

add.side.info <- function(dataset){
  dataset$side <- vector(mode="character", length = dim(dataset)[1])
  dataset$sideanum <- vector(mode="numeric", length = dim(dataset)[1])
  dataset$sidebnum <- vector(mode="numeric", length = dim(dataset)[1])
  
  
  dataset$side[dataset$avotelock_1==1] <- "A"
  dataset$side[dataset$bvotelock_1==1] <- "B"
  
  dataset$sideanum[  dataset$a33==1 |
                        dataset$a39==1 |
                        dataset$b33==1 |
                        dataset$b39==1] <-3
  dataset$sideanum[dataset$a93==1 |
                      dataset$a99==1 |
                      dataset$b93==1 |
                      dataset$b99==1] <-9
  
  dataset$sideanum[dataset$a11==1 |
                      dataset$b11==1] <-1
  
  dataset$sidebnum[  dataset$a33==1 |
                        dataset$a93==1 |
                        dataset$b33==1 |
                        dataset$b93==1] <-3
  dataset$sidebnum[dataset$a39==1 |
                      dataset$a99==1 |
                      dataset$b39==1 |
                      dataset$b99==1] <-9
  
  dataset$sidebnum[dataset$a11==1 |
                      dataset$b11==1] <-1
  
  return(dataset)
}

raw.data <- add.side.info(raw.data)

# Data for subjects that got far enought to see the
# quiz
full.data <- subset(raw.data, quizintro==1)

# Data where subjects correctly completed 
# the comprehension quiz
# This is the data used in the main manuscript
# the raw data is examined further in Supplementary Information Appendix
# to see whether the quiz affected balance
valid.data <- subset(raw.data, validation==1 )


# Simulates bargaining by taking a subset of the relative data
# pairing players/groups, tabulating the outcomes
# it is meant to be run using the bootstrap command, below.
# Depending on the "outcome" the function records different measures 
# the default is bargaining success, but it can also record
# how much members of one side earned or how accurate proposers guesses were
# about responders Minimum Acceptable Offer (MAO)
bargain.simulator <- function(thedata, indices, sideanum, sidebnum, outcome="bargaining.success"){
  # split simulation dataset into members from eligble sides based on subject's condition
  #outcome.measure <- vector(mode="numeric", length=1000)
  simdata <- thedata[indices,]
  sidea <- subset(simdata, side=="A" & sideanum==sideanum & sidebnum==sidebnum)
  sideb <- subset(simdata, side=="B" & sideanum==sideanum & sidebnum==sidebnum)
  # define where the voters are in the simulation dataset

 vb.beg <- match(x="vsideb_1",names(simdata))
 vb.end <- match(x="vsideb_21",names(simdata)) 
  # grab the simulated voters'votes
 va.votes <- sidea$vsidea
 vb.votes <- sideb[,vb.beg:vb.end]
  rr <- foreach(i = 1:1000, .combine=cbind) %dopar% { 
  
    #select/sample voters
    sidea.voters <- sample(x=1:length(va.votes), size=sideanum,replace=T)
    sideb.voters <- sample(x=1:dim(vb.votes)[1], size=sidebnum,replace=T)
    
    # tally votes
    sidea.votes <- va.votes[sidea.voters]
    
    sideb.vote.vectors <-vb.votes[sideb.voters,]
    sideb.votes <- colSums(sideb.vote.vectors, na.rm=T)
    
    # define a majority for each side
    if(sideanum>1){sidea.maj.vote <- quantile(x=sidea.votes, probs=.5)}else{
      sidea.maj.vote <- sidea.votes
    }
    if(sidebnum>1){sideb.maj.vote <- sideb.votes >= sidebnum/2+.5}else{
      sideb.maj.vote <- sideb.votes == 1
    }
    
    #find out amounts  have a majority vote
    amounts <- seq(from=0, to=40, by=2)
    
    sidea.final.offer <- sidea.maj.vote
    
    sideb.accepts <- amounts[sideb.maj.vote]
    sideb.final.offer <- sidea.final.offer%in%sideb.accepts
    
    sidea.earns <- 0
    sideb.earns <- 0
    
    if(sideb.final.offer){
      sidea.earns <- 40-sidea.final.offer
      sideb.earns <- sidea.final.offer
    }
        
    if(outcome=="bargaining.success"){return(as.numeric(sideb.final.offer))}
    if(outcome=="sideaearns"){return(c(sidea.earns))}
    if(outcome=="sidebearns"){return(c(sideb.earns))}
    if(outcome=="accuracya"){return(sqrt(mean(sidea$guessb_1 - min(sideb.accepts))^2, na.rm=T)) #rmse  
      } 
  }
  #return(mean(outcome.measure, na.rm=T))
  mean(rr)
}



# Run the simulator for groups of 5 vs 5, recording bargaining success or failure
# There was no experiment where 5 subjects actually faced another 5 subjects
# This simulation uses data from the 9 v 9 condition as part of a robustness check
# requested by the reviewers.
success55.full <-boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=5, sidebnum=5, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success55.full$t,  file="success55.full.csv")
save(success55.full, file = "success55.full.Rdata")

# Run the simulator for groups of 7 vs 7, recording bargaining success or failure
# There was no experiment where 7 subjects actually faced another 7 subjects
# This simulation uses data from the 9 v 9 condition as part of a robustness check
# requested by the reviewers.
success77.full <-boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=7, sidebnum=7, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success77.full$t,  file="success77.full.csv")
save(success77.full, file = "success77.full.Rdata")

# Run the simulator for groups of 3 vs 3, recording bargaining success or failure
success33<- boot(data=subset(valid.data, sideanum==3 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=3, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success33$t,  file="success33.csv")
save(success33, file = "success33.Rdata")

# Run the simulator for groups of 3 vs 9, recording bargaining success or failure
success39<- boot(data=subset(valid.data, sideanum==3 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=9, outcome="bargaining.success", ncpus=32, parallel="multicore")
write.csv(success39$t,  file="success39.csv")
save(success39, file = "success39.Rdata")

# Run the simulator for groups of 9 vs 3, recording bargaining success or failure
success93<- boot(data=subset(valid.data, sideanum==9 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=3, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success93$t,  file="success93.csv")
save(success93, file = "success93.Rdata")

# Run the simulator for groups of 9 vs 9, recording bargaining success or failure
success99<- boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=9, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success99$t,  file="success99.csv")
save(success99, file = "success99.Rdata")


# Run the simulator for groups of 3 vs 3, recording side A's earnings
aearns33<- boot(data=subset(valid.data, sideanum==3 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=3, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(aearns33$t,  file="aearns33.csv")
save(aearns33, file = "aearns33.Rdata")

# Run the simulator for groups of 3 vs 9, recording side A's earnings
aearns39<- boot(data=subset(valid.data, sideanum==3 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=9, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(aearns39$t,  file="aearns39.csv")
save(aearns39, file = "aearns39.Rdata")

# Run the simulator for groups of 9 vs 3, recording side A's earnings
aearns93<- boot(data=subset(valid.data, sideanum==9 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=3, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(aearns93$t,  file="aearns93.csv")
save(aearns93, file = "aearns93.Rdata")

# Run the simulator for groups of 9 vs 9, recording side A's earnings
aearns99<- boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=9, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(aearns99$t,  file="aearns99.csv")
save(aearns99, file = "aearns99.Rdata")


# Run the simulator for groups of 3 vs 3, recording side B's earnings
bearns33<- boot(data=subset(valid.data, sideanum==3 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=3, outcome="sidebearns", ncpus=4, parallel="multicore")
write.csv(bearns33$t,  file="aearns33.csv")
save(bearns33, file = "bearns33.Rdata")


# Run the simulator for groups of 3 vs 9, recording side B's earnings
bearns39<- boot(data=subset(valid.data, sideanum==3 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=9, outcome="sidebearns", ncpus=4, parallel="multicore")
write.csv(bearns39$t,  file="aearns39.csv")
save(bearns39, file = "bearns39.Rdata")

# Run the simulator for groups of 9 vs 3, recording side B's earnings
bearns93<- boot(data=subset(valid.data, sideanum==9 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=3, outcome="sidebearns", ncpus=4, parallel="multicore")
write.csv(bearns93$t,  file="bearns93.csv")
save(bearns93, file = "bearns93.Rdata")

# Run the simulator for groups of 3 vs 9, recording side B's earnings
bearns99<- boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=9, outcome="sidebearns", ncpus=4, parallel="multicore")
write.csv(bearns99$t,  file="bearns99.csv")
save(bearns99, file = "bearns99.Rdata")


# Run the simulator for groups of 3 vs 3, recording side A's accuracy at guessing side B's resevation price (RMSE)
accuracya33<- boot(data=subset(valid.data, sideanum==3 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=3, outcome="accuracya", ncpus=4, parallel="multicore")
write.csv(accuracya33$t,  file="accuracya33.csv")
save(accuracya33, file = "accuracya33.Rdata")

# Run the simulator for groups of 3 vs 9, recording side A's accuracy at guessing side B's resevation price (RMSE)
accuracya39<- boot(data=subset(valid.data, sideanum==3 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=3, sidebnum=9, outcome="accuracya", ncpus=4, parallel="multicore")
write.csv(accuracya39$t,  file="accuracya39.csv")
save(accuracya39, file = "accuracya39.Rdata")

# Run the simulator for groups of 9 vs 3, recording side A's accuracy at guessing side B's resevation price (RMSE)
accuracya93<- boot(data=subset(valid.data, sideanum==9 & sidebnum==3), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=3, outcome="accuracya", ncpus=4, parallel="multicore")
write.csv(accuracya93$t,  file="accuracya93.csv")
save(accuracya93, file = "accuracya93.Rdata")


# Run the simulator for groups of 9 vs 9, recording side A's accuracy at guessing side B's resevation price (RMSE)
accuracya99<- boot(data=subset(valid.data, sideanum==9 & sidebnum==9), statistic=bargain.simulator, R=3000, sideanum=9, sidebnum=9, outcome="accuracya", ncpus=4, parallel="multicore")
write.csv(accuracya99$t,  file="accuracya99.csv")
save(accuracya99, file = "accuracya99.Rdata")


# Load elite data from :
# LeVeck, Brad L et al. 2014. “The Role of Self-Interest in Elite Bargaining.” Proceedings of the National Academy of Sciences of the United States of America 111(52): 18536–41.

elite.data <- read.csv("pnas.elite.turk.data.csv")
elite.data <- subset(elite.data, elite==1)


# Function to simulate how well elites do at avoidng bargaining failure and how much they earn in 1 on 1 interactions.
bargain.simulator.elite <- function(thedata, indices, outcome="bargaining.success"){
  # split simulation dataset into members from eligble sides based on subject's condition
  #outcome.measure <- vector(mode="numeric", length=1000)
  simdata <- thedata[indices,]
  outcome.measure <- vector(mode="numeric", length=dim(simdata)[1])
  
  for(i in 1:dim(simdata)[1]){
    sidea <- simdata$propose[i]
    sideb <- simdata$demand[-i]
    
    if(outcome=="bargaining.success"){
      outcome.measure[i] <- mean(sidea >= sideb, na.rm=T)}
    if(outcome=="sideaearns"){
      outcome.measure[i] <- mean((100-sidea)*(sidea >= sideb), na.rm=T)*.4}
  }
  return(mean(outcome.measure, na.rm=T))
}

# Simulate elite bargaining success in Ultimatum game.
elite_success11<- boot(data=elite.data, statistic=bargain.simulator.elite, R=5000, outcome="bargaining.success", ncpus=32, parallel="multicore")
write.csv(elite_success11$t,  file="elite_success11.csv")
save(elite_success11, file = "elite_success11.Rdata")

# Simulate elite earning success success in Ultimatum game.
elite_aearns11<- boot(data=elite.data, statistic=bargain.simulator.elite, R=5000, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(elite_aearns11$t,  file="elite_aearns11.csv")
save(elite_aearns11, file = "elite_aearns11.Rdata")

## Function to simulate 1 on 1 bargaining among individuals in our experiment. Used to compared individuals to elites
bargain.simulator.oneonone <- function(thedata, indices, outcome="bargaining.success"){
  # split simulation dataset into members from eligble sides based on subject's condition
  #outcome.measure <- vector(mode="numeric", length=1000)
  simdata <- thedata[indices,]
   
  sidea <- subset(simdata, side=="A" )
  sideb <- subset(simdata, side=="B" )
  
  vb.beg <- match(x="vsideb_1",names(simdata))
  vb.end <- match(x="vsideb_21",names(simdata))
  
  outcome.measure <- vector(mode="numeric", length=dim(sidea)[1])

  for(i in 1:dim(sidea)[1]){
    va.votes <- sidea$vsidea[i]
    vb.votes <- sideb[,vb.beg:vb.end]
    if(outcome=="bargaining.success"){
      outcome.measure[i] <- mean(vb.votes[,(va.votes+2)/2], na.rm=T)}
    if(outcome=="sideaearns"){
      outcome.measure[i] <- (40-va.votes)* mean(vb.votes[,(va.votes+2)/2], na.rm=T)}
  }
  return(mean(outcome.measure, na.rm=T))
}

# Simulate bargaining success from 1 on 1 bargaining
success11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1), statistic=bargain.simulator.oneonone, R=5000, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(success11$t,  file="success11.csv")
save(success11, file = "success11.Rdata")

# Simulate earnings from 1 on 1 bargaining
aearns11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1), statistic=bargain.simulator.oneonone, R=5000,  outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(aearns11$t,  file="aearns11.csv")
save(aearns11, file = "aearns11.Rdata")

# Simulate bargaining success from 1 on 1 bargaining among experienced players (who reported having played a bargaining game before taking our experiment)
experienced_success11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1 & experience==1), statistic=bargain.simulator.oneonone, R=5000, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(experienced_success11$t,  file="experienced_success11.csv")
save(experienced_success11, file = "experienced_success11.Rdata")

# Simulate bargaining earnings from 1 on 1 bargaining among experienced players (who reported having played a bargaining game before taking our experiment)
experienced_aearns11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1 & experience==1), statistic=bargain.simulator.oneonone, R=5000, outcome="sideaearns", ncpus=32, parallel="multicore")
write.csv(experienced_aearns11$t,  file="experienced_aearns11.csv")
save(experienced_aearns11, file = "experienced_aearns11.Rdata")

# Simulate bargaining success from 1 on 1 bargaining among inexperienced players (who reported they had not played a bargaining game before taking our experiment)

inexperienced_success11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1 & experience!=1), statistic=bargain.simulator.oneonone, R=5000, outcome="bargaining.success", ncpus=4, parallel="multicore")
write.csv(inexperienced_success11$t,  file="inexperienced_success11.csv")
save(inexperienced_success11, file = "inexperienced_success11.Rdata")

# Simulate bargaining earnings from 1 on 1 bargaining among inexperienced players (who reported they had not played a bargaining game before taking our experiment)
inexperienced_aearns11<- boot(data=subset(valid.data, sideanum==1 & sidebnum==1 & experience!=1), statistic=bargain.simulator.oneonone, R=5000, outcome="sideaearns", ncpus=4, parallel="multicore")
write.csv(inexperienced_aearns11$t,  file="inexperienced_aearns11.csv")
save(inexperienced_aearns11, file = "inexperienced_aearns11.Rdata")



#############################
# Code for Figures          #
#############################

# Load recorded values from prior
# simulations and make labels

bargaining.success <- read.csv("sim_results.bargaining.success.csv")
bargaining.success$totalnum <- c(2,2,2,6,12,12,18)
bargaining.success$Condition <- factor(bargaining.success$Condition, 
                                       levels = c("c11_inexperienced",
                                                  "c11_experienced",
                                                  "c11_elite",
                                                  "c33",
                                                  "c39",
                                                  "c93",
                                                  "c99"),
                                       labels=c("11\ninexperienced\nsubjects",
                                                "11\nexperienced\nsubjects", 
                                                "11\ninternational\npolicy elite",
                                                "33",
                                                "39",
                                                "93",
                                                "99"))

bargaining.success.main <- bargaining.success[-1:-3,]

# Figure 1a
bargaining.success.barplot <- ggplot(bargaining.success.main, aes(y=Bargaining.Success, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=Bargaining.Success-SE, ymax=Bargaining.Success+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + scale_y_continuous(expand = c(0, 0), limits=c(0,1.1))+ylab("% Bargaining Successes") +xlab("Experimental Condition") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("Autocracy (P)\nAutocracy (R)", "Autocracy (P)\nDemocracy (R)", "Democracy (P)\nAutocracy (R)", "Democracy (P)\nDemocracy (R)"))
bargaining.success.barplot
ggsave(bargaining.success.barplot, file="bargaining.success.barplot.pdf", width=8, height =6)

bargaining.success.experts <- bargaining.success[-4:-6,]

# Figure 2a
bargaining.success.experts.barplot <- ggplot(bargaining.success.experts, aes(y=Bargaining.Success, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=Bargaining.Success-SE, ymax=Bargaining.Success+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + scale_y_continuous(expand = c(0, 0), limits=c(0,1.1))+ylab("% Bargaining Successes") +xlab("Experimental Condition") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("Inexperienced\nIndividuals", "Experienced\nIndividuals", "International\nPolicy Elites", "Democratic\nDyads"))
bargaining.success.experts.barplot



earnings<- read.csv("sim_results.earnings.csv")
earnings$totalnum <- c(2,2,2,6,12,12,18)
earnings$Condition <- factor(earnings$Condition, 
                                       levels = c("c11_inexperienced",
                                                  "c11_experienced",
                                                  "c11_elite",
                                                  "c33",
                                                  "c39",
                                                  "c93",
                                                  "c99"),
                                       labels=c("11\ninexperienced\nsubjects",
                                                "11\nexperienced\nsubjects", 
                                                "11\ninternational\npolicy elite",
                                                "33",
                                                "39",
                                                "93",
                                                "99"))

earnings.main <- earnings[-1:-3,]

# Figure 1b
earnings.barplot <- ggplot(earnings.main, aes(y=Earns, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=Earns-SE, ymax=Earns+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + scale_y_continuous(expand = c(0, 0), limits=c(0,20))+ylab("Mean Earnings (cents)") +xlab("Experimental Condition") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("Autocracy (P)\nAutocracy (R)", "Autocracy (P)\nDemocracy (R)", "Democracy (P)\nAutocracy (R)", "Democracy (P)\nDemocracy (R)"))
earnings.barplot
ggsave(earnings.barplot, file="earnings.barplot.pdf", width=8, height =6)


earnings.experts <- earnings[-4:-6,]

# Figure 2b
earnings.experts.barplot <- ggplot(earnings.experts, aes(y=Earns, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=Earns-SE, ymax=Earns+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + scale_y_continuous(expand = c(0, 0), limits=c(0,20))+ylab("Mean Earnings (cents)") +xlab("Experimental Condition") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("Inexperienced\nIndividuals", "Experienced\nIndividuals", "International\nPolicy Elites", "Democratic\nDyads"))
earnings.experts.barplot


# Plot belief about MAO vs Proposers Vote Figure A2 panel A in SI Appendix
beliefa.plot <- ggplot(data=subset(valid.data, side=="A"), aes(x=guessb_1, y=vsidea))+geom_jitter(aes(), alpha=.3)+geom_smooth(, method="lm")+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + xlab("Belief About Responders' \nMinimum Acceptable Offer")+ ylab("Proposers' Vote on\n Amount to Offer")+xlim(c(0,40))+ylim(c(0,40))
beliefa.plot <- ggplot(data=subset(valid.data, side=="A"), aes(x=guessb_1, y=vsidea))+geom_jitter(aes(), alpha=.3)+geom_smooth(, method="loess", span=1)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + xlab("Belief About Responders' \nMinimum Acceptable Offer")+ ylab("Proposers' Vote on\n Amount to Offer")+xlim(c(0,20))+ylim(c(0,40))
beliefa.plot


# Plot own demand as a responder vs offer as a proposer  
## Figure A2 panel B in SI Appendix
#using data from:
# LeVeck, Brad L et al. 2014. “The Role of Self-Interest in Elite Bargaining.” Proceedings of the National Academy of Sciences of the United States of America 111(52): 18536–41.

undergrads <- read.csv("elite.ug.csv")

undergrads <- subset(undergrads, elite==0)

offer.demand.plot <- ggplot(data=undergrads, aes(x=ult.min*.4, y=ult.propose*.4))+geom_jitter(aes(), alpha=.3)+ geom_abline(intercept = 0, slope=1, color="grey")+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + xlab("Proposers' Demand When Playing\n as the Responder in a Separate Game")+ ylab("Proposers' Offer")+xlim(c(0,40))+ylim(c(0,40))
offer.demand.plot 


# Create a bar plot on how accurate subjects 
# were at guessing the MAO of side B

belief.accuracy<- read.csv("belief.accuracy.csv")
belief.accuracy$totalnum <- c(6,12,12,18)
belief.accuracy$Condition <- factor(belief.accuracy$Condition, 
                             levels = c("c33",
                                        "c39",
                                        "c93",
                                        "c99"),
                             labels=c("33",
                                      "39",
                                      "93",
                                      "99"))

belief.accuracy.barplot <- ggplot(belief.accuracy, aes(y=RMSE, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=RMSE-SE, ymax=RMSE+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) +ylab("RMSE of Proposer's Belief") +xlab("Experimental Condition") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("Autocracy (P)\nAutocracy (R)", "Autocracy (P)\nDemocracy (R)", "Democracy (P)\nAutocracy (R)", "Democracy (P)\nDemocracy (R)"))
belief.accuracy.barplot
ggsave(belief.accuracy.barplot, file="belief.accuracy.barplot.pdf", width=8, height =6)



# Analysis to see if offers were fundamentally different across conditions. 

median(valid.data$vsidea[valid.data$side=="A" & valid.data$sideanum==3], na.rm=T)
mean(valid.data$vsidea[valid.data$side=="A" & valid.data$sideanum==3], na.rm=T)

median(valid.data$vsidea[valid.data$side=="A" & valid.data$sideanum==9], na.rm=T)
mean(valid.data$vsidea[valid.data$side=="A" & valid.data$sideanum==9], na.rm=T)


# Examine attrition from the quiz
# as reported in the Supplementary Information Appendix
# section titled Recruitment of Subjects and Analysis of Balance
raw.data$qnumA <-raw.data$qnumA3_1
raw.data$qnumA[is.na(raw.data$qnumA)] <- raw.data$qnumA9_1[is.na(raw.data$qnumA)]
raw.data$qnumA[is.na(raw.data$qnumA)] <- raw.data$qnumA1_1[is.na(raw.data$qnumA)]

raw.data$qnumB <-raw.data$qnumB3_1
raw.data$qnumB[is.na(raw.data$qnumB)] <- raw.data$qnumB9_1[is.na(raw.data$qnumB)]
raw.data$qnumB[is.na(raw.data$qnumB)] <- raw.data$qnumB1_1[is.na(raw.data$qnumB)]

raw.data.test <- subset(raw.data, quizintro==1 & !(is.na(vsidea) & side=="A"))
raw.data.test$validation[is.na(raw.data.test$validation)]<-0
raw.data.test$validation2 <- raw.data.test$validation


raw.data.test$validation2[!is.na(raw.data.test$qnumA) & raw.data.test$qnumA != raw.data.test$sideanum] <-1
raw.data.test$validation2[!is.na(raw.data.test$qnumB) & raw.data.test$qnumB != raw.data.test$sidebnum] <-1


raw.data.test$condition <- paste(raw.data.test$side,raw.data.test$sideanum, raw.data.test$sidebnum,sep="")

# Check if those who passed the quiz are different from theoretical distribution
# As shown in Supplementary Information Appendix Table A1
chisq.test(table(raw.data.test$condition[raw.data.test$validation==1]), p = c(1/14, 1/14, 1/14, 3/14, 1/14, 1/14, 1/14, 3/14, 1/14, 1/14))

# Check if those who passed the quiz are different from theoretical distribution
# as long as they remembered the number of people on the other side
# As shown in Supplementary Information Appendix Table A2
chisq.test(table(raw.data.test$condition[raw.data.test$validation2==1]), p = c(1/14, 1/14, 1/14, 3/14, 1/14, 1/14, 1/14, 3/14, 1/14, 1/14))

#############################
# Additional analyses for R&R
#############################
library(dplyr)
library(ggplot2)

# read in Polity IV data

polity4 <- read.csv("polity4.csv", stringsAsFactors = F)

polcon <- read.csv("polcon2012.csv", stringsAsFactors = F) %>%
  select(cyear, polconiii)

polity4 <- polity4 %>%
  mutate(regime = ifelse(polity2 > 5,"Democracy",
                            ifelse(polity2 < 6 & polity2 > -6, "Anocracy",
                                     ifelse(polity2 < -5, "Autocracy", NA)))) %>%
  mutate(regime = factor(regime, levels = c("Democracy", "Anocracy", "Autocracy") ))
pol_merged <-  left_join(x=polity4, polcon, by = "cyear")


polity_summary <- pol_merged %>%
  group_by(regime) %>%
  summarise(mean_xconst = round(mean(xconst, na.rm =T), 2),
            mean_parcomp = round(mean(parcomp, na.rm=T), 2),
            mean_polconiii = round(mean(polconiii, na.rm=T), 2)
            ) %>%
  na.omit()
polity_summary
polity_summary


# Graph from SI Appendix A1
interpolated_bargaining.success <- read.csv(file="interpoltated_sims.bargaining.success.csv")
interpolated_bargaining.success$totalnum <- c(6,10,14,18)
interpolated_bargaining.success$Condition <- factor(interpolated_bargaining.success$Condition, 
                                       levels = c("c33",
                                                  "c55",
                                                  "c77",
                                                  "c99"),
                                       labels=c("33",
                                                "55",
                                                "77",
                                                "99"))


interpolated_bargaining.success.barplot <- ggplot(interpolated_bargaining.success, aes(y=Bargaining.Success, x=Condition))+geom_bar(aes(fill=totalnum), stat="identity", color=I("black"), show_guide=F)+geom_errorbar(aes(ymin=Bargaining.Success-SE, ymax=Bargaining.Success+SE, width=.06),  size=.7, show_guide=F)+theme(axis.title.x=element_text(size=20, vjust=-.5, face="bold"), axis.title.y=element_text(size=20,angle=90, vjust=0.8, face="bold"), axis.text.y=element_text(size=16,angle=90, vjust=-.5, hjust=0.55, face="bold") , axis.text.x=element_text(size=16, face="bold"), panel.grid.major=element_blank(), panel.grid.minor=element_blank(), panel.background=element_blank(), plot.margin=unit(c(1,1,1,1),"cm"), axis.line=element_line(), legend.position=c(0,1), legend.justification=c(0,1), legend.title= element_text(size=14, face="bold"), legend.text=element_text(size=12, face="bold")) + scale_y_continuous(expand = c(0, 0), limits=c(0,1.1))+ylab("% Bargaining Successes") +xlab("Group Size") + scale_fill_gradient2(low='white', high='red')+guides(fill=FALSE)+scale_x_discrete(labels=c("3 (P)\n3 (R)", "5 (P)\n5 (R)", "7 (P)\n7 (R)", "9 (P)\n9 (R)"))
interpolated_bargaining.success.barplot
ggsave(interpolated_bargaining.success.barplot, file="interpolated_bargaining.success.barplot.pdf", width=8, height =6)

# Check that democracies are not earnign substantially less
# conditional on bargaining success
# assumes you have already run and saved the appropriate simulations
# Conditional earnings
# all are within less than 1 cent of one another
load("cond_aearns99.Rdata")
load("cond_aearns93.Rdata")
load("cond_aearns39.Rdata")
load("cond_aearns33.Rdata")

mean(cond_aearns33$t,na.rm=T)
sd(cond_aearns33$t,na.rm=T)
mean(cond_aearns39$t,na.rm=T)
sd(cond_aearns39$t,na.rm=T)
mean(cond_aearns93$t,na.rm=T)
sd(cond_aearns93$t,na.rm=T)
mean(cond_aearns99$t,na.rm=T)
sd(cond_aearns99$t,na.rm=T)